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ABSTRACT 

A general method to find an effective potential of interaction between far 
separated 2D and 3D solitons is elaborated, including the case of 2D vortex 
solitons. The method is based on explicit calculation of the overlapping term 
in the full Hamiltonian of the system (without assuming that the "tail" of 
each soliton is not affected by its interaction with the other soliton). The 
result is obtained in an explicit form that does not contain an artificially 
introduced radius of the overlapping region. The potential applies to spatial 
and spatiotemporal solitons in nonlinear optics, where it helps to solve various 
dynamical problems: collisions, formation of bound states (BS's), etc. In 
particular, an orbiting BS of two solitons is always unstable. In the presence 
of weak dissipation and gain, the effective potential can also be derived, 
giving rise to bound states similar to those recently studied in ID models. 

PACS numbers: 42.50.Rh; 42.50.Ne; 42.65.Vh; 52.35.Sb; 03.40.Kf 
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1 INTRODUCTION 

Recent progress in studies of two-dimensional (2D) solitons in models of 
non-Kerr nonlinear optical media has attracted a lot of interest to their in- 
teractions. 2D vortex solitons and interactions between them in the quintic 
nonlinear Schrodinger (QNLS) equation were studied in JL], nonplanar inter- 
actions between 2D solitons in a medium with the quadratic (x^) nonlinear- 
ity were considered, numerically and analytically, in , and various features 
of the interaction between 2D solitons in photorefractive media were revealed 
by numerical simulations and direct experiments p, ^, |5|]. The nonlinearity 
must be non-Kerr because the usual cubic (Kerr) self-focusing term gives rise 
to collapse in 2D and 3D (three-dimensional) cases. As it was demonstrated 
in ||, the collapse does not take place in any physical dimension in the model 
with the x nonlinearity. This opens way to stable 2D and 3D spatiotempo- 
ral solitons, or "light bullets" (LBs) 0. The x (2) LBs were recently studied 
in detail in || and |§. 

The objective of this work is to find an effective potential of interac- 
tion between 2D and 3D solitons in isotropic media (note that, as it was 
demonstrated in a very recent experimental work M, the interaction of 2D 
solitons in intrinsically anisotropic photorefractive media is, in effect, prac- 
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tically isotropic too). The interaction potential is necessary to solve various 
dynamical problems, such as collisions, formation of bound states of soli- 
tons, etc., including a practically important problem of designing all-optical 
switching by means of interaction between 2D optical solitons. It will be 
demonstrated that a universal effective potential can be obtained analyti- 
cally by means of a technique which generalizes that developed for the ID 



solitons in ||10|| . As a paradigm model, one can take the multidimensional 
quintic Ginzburg- Landau (QGL) equation, 

iv t H — V 2 f + \v \ 2 v — a\v |S = —iv + ryiV 2 t> + Z7 2 |f | 2 f — 173 |i> | 4 f, (1) 

where the coefficients a and 71,2,3 are positive. The QNLS equation is a con- 
servative version of flU), without its right-hand side. The quintic defocusing 
term ~ a is included in order to prevent the collapse. Note that this term 
is not merely the simplest one that stabilizes the model: according to ex- 



perimental data [IT], the combination of the focusing cubic and defocusing 
quintic terms adequately models the nonlinear optical properties of some real 
materials. The first two terms on the right-hand side of (|TJ) take into regard 
linear losses, the cubic term 72 accounts for nonlinear gain which compen- 
sates the losses, and the quintic dissipation term ~ 73 provides for the overall 



stabilization of the model. The QGL equation was first introduced in [12 



(in the 2D form), and its ID (one-dimensional) version has later attracted 



a great deal of interest (see, e.g., Ej] and references therein). In particular, 



stable localized pulses in the ID QGL equation were found in [14]] for the 
case of weak dissipation (relevant for the applications to nonlinear optics), 
< 71,2,3 "C 1, that will be also assumed here. The existence of the stable 
pulses in the opposite limit of strong dissipation was independently shown 
in three different works [|ll|. Actually, the model (H) is selected just for the 
reference, as the one that certainly gives rise to stable multidimensional soli- 
tons; as it will be seen below, the derivation of the effective potential for the 
interaction between the solitons, presented in this work, is quite universal 
and may be applied to any conservative or weakly dissipative model that 
supports multidimensional solitons. 

Note that stable 2D solitons, as well as two-soliton bound states, were also 
found numerically in a model with the quintic nonlinearity similar to that in 
Eq. ([!]), in which, however, the linear part is of a higher order, containing 
the operators d 2 /dt 2 and V 4 |16[]. However, that model is essentially more 
complicated than (|l]), and its physical applications are less clear. 

The paper is organized as follows. In section 2, the 2D and 3D soliton so- 
lutions are briefly considered, with emphasis on the form of their asymptotic 
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"tails", which determine the effective interaction potential. In the same sec- 
tion, the model ([[]) is also reformulated in terms of nonlinear optics, where it 
finds applications of two types: the description of spatial cylindrical solitons 
in the bulk medium, and 2D and 3D LBs in the, respectively, 2D nonlin- 
ear waveguide or 3D bulk. The multidimensional solitons in the model with 
the x nonlinearity have their own peculiarities, which are summarized in 
a separate subsection in section 2. In section 3, the interaction potential is 
analytically derived, in a general form, for the 2D and 3D solitons. In the 
same section, the interaction potential for LBs is also considered. In par- 
ticular, the potential may be spatiotemporally anisotropic for the x LBs, 
while in the other models it can always be cast into an effectively isotropic 
form. Concluding remarks are collected in section 4, including a discussion 
of a possibility of existence of bound states of the solitons. In particular, it 
is concluded, in accord with the recent results obtained for the x^ spatial 
solitons in || , that a bound state of two solitons orbiting around each other 
may exist (in the dissipationless model), but it is always unstable. In the 
presence of the weak dissipation and gain, there are bound states of quiescent 
solitons, quite similar to those recently studied in the ID model (that may 
be both unstable and almost stable). New possible states in the 2D and 3D 
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cases are soliton lattices and "molecules" . 



2 TWO- AND THREE-DIMENSIONAL SOLI- 
TONS 

2.1 The General Case 

A general stationary solution to Eq. ([[]) is v — exp (—iuit) V^(r), where ^(r) 
satisfies the equation 

-V 2 V+\V\ 2 V-a\V\ 4 V + uV = -iV + z 7l VV + i^ 2 \V\ 2 V -i^V^V. (2) 

In the 2D case, the solution is restricted to the form 

V(r) = exp (is9) V(r), s = 0, ±1, ±2, ... , (3) 

where r and 6 are the polar coordinates, s ^ corresponding to a vortex 
soliton, and V(r) exponentially decays at r — > oo. From the consideration of 
Eq. ([I]) it follows that the asymptotic form of the soliton at r — > oo is 

V(r) A,r~ 1/2 exp (-/cr) , (4) 



and, at r — > 0, 

V(r) w a s r |s| (6) 

(i.e., the vortex soliton has a hole in its center), with unknown constants 
A s and a s . The expansion of k in ([|) employs the fact that, in the weakly 
dissipative regime, 71 is small, and to ^> 1, as the dissipation coefficient in 
front of the term — v in Eq. (|TJ) is 1 . 

The stability of the s = soliton in the model (fl|) is very plausible, and, 
in the conservative version of (H), the stability of the vortex soliton with 
\s\ = 1 was numerically demonstrated in 0. It is not known if the solitons 
with \s\ > 1 are stable (note that all the bright vortex solitons are unstable 
in the model, see, e.g., [17]). Below, an arbitrary integer value of s will 



be kept, as the potential can be derived in the general case, provided that 
the two solitons have s± — ±S2- 

Description of 3D solitons with the internal "spin" is a rather complicated 
problem, therefore only the 3D solitons with the zero spin will be considered 
here. The corresponding solution is sought for in the form of Eq. (^) with 
s — 0, and with the difference that r is now the radial variable in the 3D 
space, hence 

V(r) ps A s r~ l exp (— nr) (7) 



at r — > oo. 

In the conservative version of the model, the frequency uj < is an 
arbitrary parameter of the soliton, while the amplitude A s , that can be found 
numerically, is a function of uj (as well as a s in Eq. @). In the presence 
of the weak dissipation and gain, an actual soliton solution is selected from 
the continuous family as the one providing for a balance of the "number of 
photons", Jq 00 \ V{j-)\ 2 r D ~ l dr [14]. In that case, the value of uj should also be 
found numerically. Below, uj and A s will be treated as given parameters. 

In the application to the nonlinear optics, the QNLS version of the 2D 
model ([!]) describes time-independent light distributions in a 3D medium, so 
that the variable t is not time, but the propagation coordinate. The dynamics 
of LBs in 2D and 3D optical media is governed by an equation that is also 
similar to Eq. (|l|). Neglecting the dissipative part, the corresponding QNLS 
equation is 

iv x + \ (Vju + v TT ) + \v\ 2 v - a\v\ 4 v = 0, (8) 
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where v is the envelope of the electromagnetic waves, z and r = t — z/c gr 
are the propagation coordinate and the so-called retarded time, c gr being 
the mean group velocity of the carrier wave, and the operator acts on 
the transverse coordinate(s). In Eq. (||), anomalous temporal dispersion 
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(accounted for by the term v TT ) is assumed. A spatiotemporal-soliton solution 
to Eq. dS|) (i.e., LB) can be sought for in the form (cf. Eq. (0)) 

v = exp (ikz) V(0, f = + r 2 . (9) 

Here, k is the propagation constant and rj_ is the transverse coordinate in 
the 2D model, or the radial variable in the transverse plane in the 3D model. 
In the latter case, a more general solution with a "spatiotemporal spin" can 
be looked for in the form 

v = exp (ikz + isd) V(f), (10) 

where this time 9 is the formal angular coordinate on the plane (r±,r). The 
solution ([[(]) has a "hole" in its center, cf. Eq. @. The asymptotic form of 
all the LB solutions at £ — >• oo is similar to that given above by Eqs. (^), (|5D 
and (0): 

V(O^Ar (D - 1)/2 exp(-^). (11) 

2.2 The model with the quadratic nonlinearity 

An allied physically important model is that describing multidimensional 
media [§]: 

iV z + g ( V ± W + V rr) - V + V*W = 0, (12) 
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2iw z + - (V 2 ± w + 5w TT ) - + \ 2 = 0, (13) 

where v and w are envelopes of the fields at the fundamental harmonic (FH) 
and second harmonic (SH), 7 > is a mismatch parameter, and 5 is a relative 
coefficient of the temporal dispersion. In the real physical situations, 5 < 1, 
including negative values (which correspond to the normal dispersion at the 
second harmonic). As it was shown in || and ||, the spatiotemporal soliton 
solutions to Eqs. ( |T2| ) and fll3|) can only exist if S > 0. However, the solution 
cannot be sought for in the form (|9]), except for the unrealistic special case 
5 = 1. The asymptotic form of the soliton can be nevertheless easily found 
from the linearized versions of Eqs. ( |T2D and flT3"|), cf. Eqs. (|])-(0) and (Q): 



Ar (D - 1)/2 exp(-v / 2e), (14) 
5r (D - 1)/2 exp f-J^t) , ? = r\ + rV. (15) 



Here, only the case s = is considered, and the propagation constant is not 
explicitly introduced, as it may be absorbed by the mismatch parameter 7. 

The consideration of Eq. (|T3|) readily demonstrates that, while the asymp- 
totic expression ( fH] ) for FH is always relevant, the expression fll5"D makes 
sense only if it decays at r±,r — > 00 not faster than v 2 . Further straight- 
forward analysis shows that this condition is always satisfied, provided that 
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7 < 45, and never satisfied, if 7 > 4. In the intermediate case 

45 < 7 < 4 (16) 

(recall that the physical constraint is 5 < 1, and 7 = 4 has a special meaning 
corresponding to the exact matching between FH and SH ||), the condition 
holds, on the plane of the variables r and r, inside the sector 

{r/r^f < i^U (17) 
7 — 4o 

and does not hold outside this sector. The solitons for which this con- 
dition holds, i.e., the shape of their asymptotic "tails" in each harmonic is 
independently determined by the corresponding linearized equations, may be 
naturally called free-tail solitons. 

In the case when the above condition does not hold, i.e., as given by 



Eq. flUf) , decays slower than w s as per Eq. (JLaj ), the latter expression does 
not apply. In this case, the actual decay of the SH tail is governed by the 
quadratic term in Eq. (13), while Eq. ( 14j ) remains valid. The final result 



which, in this case, replaces Eq. (|T^) is 

(note that the underlying conditions that define the present case guarantee 
that the expression flTlf ) is not singular). The solitons of this type may be 
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called, on the contrary to the free-tail ones defined above, tail-locked solitons. 
Note that, in the intermediate case (|16|), the free-tail asymptotic expression 
dUD holds inside the sector fll7f) , while the tail- locked expression (|18|) is valid 
outside the sector. 



3 THE INTERACTION POTENTIAL 

3.1 The General Case 

Coming back to the paradigm model (|1|), one notes that the conservative (left- 
hand) side of its stationary version (0) can be derived from the Hamiltonian 

H = J(l IWI 2 - \\V\ A + la\V\ 6 - u;|V| 2 ) dr. (19) 
The Hamiltonian allows one to define an effective interaction potential for two 



separated solitons [18|]. In the original works, the wave field corresponding 
to the two-soliton configuration was postulated to be a linear superposition 
of two isolated solitons. This was substituted into the Hamiltonian, and a 
term produced by the overlapping of the "body" of each soliton with a weak 
"tail" of the other one was identified as an effective interaction potential. 
This approach requires actual calculation of the corresponding integral term 
in flT9|), a drawback being that a distortion of the "tail" due to its interaction 
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with the other soliton is ignored. In this work, a more consistent approach 
will be developed, following that elaborated for the ID solitons in jnj. The 
method is based on representing the wave field in a vicinity of each soliton 
in the form 

v{r,t) = exp {-iut) %{r) + V t {r)\ , (20) 

where V s (r) is the isolated soliton (|3|), V t (r) is a small tail generated by the 
second soliton, and the influence of a given soliton on the other soliton's "tail" 
is not neglected. The distance R between the centers of the two solitons is 
assumed to be essentially larger than the soliton's size ~ re -1 , see Eq. (f|). A 
similar structure of the wave field is assumed near the center of the second 
soliton. 

Only the case when the interacting solitons are identical is considered 
(in particular, s\ = ±S2 for the 2D vortex solitons, and the amplitudes \A S \, 
defined as per Eq. (f|), are equal), hence the solitons have the same frequency 
lu, which allows one to define a phase difference ip between them. The case of 
the identical solitons is the most relevant one, as parameters of the solitons 
are, in a real physical situation, uniqely selected by the above-mentioned 
balance between the gain and dissipation. 



The next step is, as it was said above, to insert the expression (p0|) into 
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the Hamiltonian fll^D and calculate the overlap term in an area around the 
first soliton, adding then a symmetric contribution from the vicinity of the 
other soliton. In the first approximation, only the terms linear in V t are to 
be kept, which yields the following expression for the effective interaction 
potential, 

u d (r,i/>) = [J i^vv s -vv;-\v s \ 2 v s v; + a\v s \ 4 v s v;-uv s v t ^dr 

+c.c.] + {l # 2}, (21) 

where the subscript D pertains to the dimension. Here, c.c. stands for the 
complex conjugate expression, the integration is assumed over the overlap- 
ping region in a vicinity of the first soliton, and {1 # 2} is the symmetric 
contribution from the second soliton. Applying the Gauss theorem to the 
first term in (|2l|) , one transforms, in the 2D case, the expression into 
the form 

U D (R,i>) = {[- J (^V 2 V s + \V s \ 2 V s -a\V s \*V s + ujV^V t *dr 
+^l/;(n-V)P]+c.c.} + {l # 2}, (22) 

where the surface integral term is taken over a closed contour surrounding 
the first soliton, n being a local vector normal to the contour. As the contour, 
one can choose a circumference whose center coincides with that of the first 
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soliton (Fig. 1). The radius p is chosen so that 

K~ l < p < (23) 

i.e., it is much larger than the size of the soliton, but much smaller than the 
separation between the two solitons. The final objective will be to obtain an 
expression that does not depend on the auxiliary radius p. To this end, it 
will be necessary to supplement the condition (|23| ) by the additional one 

p 2 /R < k~\ (24) 

which is obviously compatible with ([2"3"|). 

In the 3D case, the difference is that the surface integral in Eq. ( P?2"D is 
taken over a sphere of the radius p, so that Fig. 1 shows the central cross 
section of the 3D situation. The conditions (|23|) and ( ^4|) pertain equally 
well to the 3D case. 

At this stage of the analysis, the dissipative terms in Eq. (Q) are still 
neglected. Because V s is an exact single-soliton solution to Eq. (0), the 
first integral term in Eq. ( |22"D vanishes. The conditions (|23| ) allow one to 
substitute both V t and V s in the surface integral term in Eq. (^) by the 
asymptotic expressions (|]), which yields, in the 2D case, 

U 2 {R,ip) = - ] J^\A S \ 2 y/pe~ Kp [J^ r~ 1/2 exp (iip +is x - is 2 v) exp (-nr) d6 
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+c.c.] + {l # 2}. (25) 
The angles rj and 9 and the radius 



r = y/{R + p cos 6) 2 + p 2 sin 2 6 = R + pcos6 + ^ (p 2 /i?) sin 2 + ... (26) 

are defined in Fig.l, the condition (p3|) being used to expand the radical 
in ( p6|) . Substituting the expansion into fl25|) and taking into regard the 
condition (|24|), in the first approximation it is enough to keep the first two 
terms from (|26| ) in exp (— nr), and only the first term in r -1 / 2 . Additionally, 
in the same approximation one may set rj = 0, which leads to 

U 2 (R,il>) = - ] P^\A S \ 2 R~ 1/2 ^e~ Kp [e- KR J** exp (iip + istf) exp (-k P cos 6) d9 
+c.c.] + {l # 2}. (27) 

The integral in ( p7|) can be calculated exactly in terms of the Bessel 
functions, but this is not necessary. Indeed, taking into regard, in line with 
the previous approximations, that Kp 3> 1, the Laplace approximation can 
be applied to the integral, a dominant contribution coming from a vicinity 
of the point 9 = n (point A in Fig. 1): 

exp (—Kp cos 6) d9 w v / 2^(/tp) _1/2 e +Kp . (28) 

Jo 



Substituting this into (|27|) , one sees that the p-dependent multipliers in (g7j) 
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are exactly cancelled by p~ 1 / 2 e +Kp from (|2"8|). This cancellation (in the lowest- 
order approximation considered here) is a crucial result, as it makes the 
effective potential independent of the auxiliary radius p. 

Of course, the dependence on p will not disappear if one tries to calcu- 
late higher-order corrections (with respect to -R -1 ) to the effective potential. 
Actually, this implies that the effective interaction potential, treating the 
solitons as particles, can be consistently defined only in the lowest-order ap- 
proximation. At the higher orders, it is necessary to explicitly take into 
regard deformation of the solitons by the interactions, which is not an objec- 
tive of the present work. 

In the term {1 ^ 2} in fl2"T|), ip is replaced, according to its definition, by 
—ip, and the dominant point in the surface integral is 9 = 0. This means that 
the term {1 ^ 2} is obtained by the change ip —>■ —ip, siti — > —sin. Finally, 
in the multiplier e~ KR in (^7|), small q = — Iran should be also taken into 
regard (see Eq. (||)), as it gives rise to an important effect, viz., long-period 
oscillations in the exponentially decaying tail of the interaction potential (10 



Note that the potential does not directly take into account the model's small 
dissipative part; however, that part indirectly affects the potential, inducing 
the oscillations in the solitons' tails in (01) through Imre. 
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With regard to what was said above, the final expression for the potential 

1 /2 

U 2 (R,tp) = -2v / 2^|A s | 2 (-l) s cosV (y^/R) exp(- v /z 2^ J R)cos(g J R), 

(29) 

where s is either s± or S2 = ±si, both giving the same value. Except for 
the factors (\/—2(j/R\ ^ and (— l) s , which are specific for the 2D case, 
the potential (|29|) is essentially the same as that obtained in the similar ID 



models in flO|. 



In the 3D case, the consideration is also limited to the interaction of iden- 
tical solitons (as it was said above, only the spinless solitons are considered 
in the 3D case). The above expression fl2"2] ) yields the interaction potential in 
the 3D case too (recall that, in this case, the integration in the surface term 
is over the sphere). As well as in the 2D case, the first term in (|22]) vanishes 
in the approximation that neglects the direct influence of the dissipation, 
and the integration over the sphere is dominated by a contribution from a 
small vicinity of the point A (Fig. 1). Substituting into Eq. (|22|) the 3D 
asymptotic expressions (0) for V s and V t and the expansion (p6|), one arrives, 
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instead of the integral (P5|), at 
2?r jT exp (-Kpcos9) sine d6 = 2ti(k P )- 1 (e +Kp - e" Kp ) « 27r(>p)- 1 e 



-1„+KP 



(30) 

With regard to (|30|), the final expression for the effective interaction potential 



in the 3D case becomes (cf. (p9|) ) 

U 3 (R,i()) = -4vr|A,| 2 cos^ R^expl-V^jR) cos(qR). (31) 



Note that the auxiliary radius p is completely canceled out in the final ex- 
pressions ( p9|) and (|31|). 

The potentials (|29| ) and (^1|) can be as well applied to the description of 
the interaction between the 2D and 3D spatiotemporal solitons (LBs), given 
by Eqs. (H) and (£|), (0). The differences from the above results are that 
q = (recall the dissipation was completely neglected in the LB models), 
uj must be replaced by —k, and the separation R between the solitons is 
replaced by the spatiotemporal separation E defined according to Eq. (|): 



R\ + T\ (32) 

R± and T being, respectively, the separation between the solitons in the 
transverse direction and the temporal delay between them. 
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3.2 The Model with the Quadratic Nonlinearity 

The interaction potential for the x^ solitons has its own peculiarities. For 
the spatial (stationary) 2D x^ solitons, the exponentially decaying poten- 
tial with two components, generated by the FH and SH fields, was postu- 
lated in 0. The interaction between the x LBs is more complicated, be- 
cause the nonstationary model (|1^) and ( |T3| ) is, effectively, spatiotemporally 
anisotropic, as it was explained in detailed in the previous section, see Eqs. 
(T|), flToD, and (|T8|). A straightforward consideration demonstrates that, in 



both 2D and 3D cases, the SH-generated interaction potential dominates at 
7 < 5, so that the potential is given by Eqs. (|29|) and (|3lD , with u replaced 



by —7, q = 0, and R replaced by E = J R± + 5 _1 T 2 , cf. Eqs. (p~5f ) and 



([32]). On the contrary to this, at 7 > 1, the FH-generated interaction always 



dominates, which means that one should use the potentials ( P9|) and (pl|), 
with uj = — 1 and R replaced by H defined as per Eq. d32[) . In the interme- 
diate case 5 < 7 < 1 (cf. Eq. (|T6|); recall that the physically relevant case 
is 5 < 1), the interaction potential turns out to be truly anisotropic in the 
plane (R±,T): the SH-generated interaction dominates inside the sector (cf. 
Eq. (0)) 



(T/ Rl _Y < L^Lj, (33) 
7 — 
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and the FH-generated interaction dominates outside the sector (^). Ac- 
cordingly, one should substitute R in the expressions fl29|) and for the 
interaction potential by 5 inside the sector (|33|) , and by H outside of it. 

4 CONCLUDING REMARKS 

The effective interaction potentials (|29|) and ( pi]) can give rise to bound states 
(BS's) of two solitons. In the presence of the dissipation and gain, it makes 
sense to consider only BS's of quiescent solitons, as any motion is suppressed 
by a friction force. Because the form of the potentials is essentially the same 
as in ID, the situation is not different from the ID case, which was recently 
studied in detail [JTJJ. There are two types of BS's, with the phase difference 
between the solitons ip = or ir, and with ip = tt/2. The BS's of the former 
type are saddles, while the BS's of the latter type have imaginary stability 
eigenvalues. The fact that the BS's with ip = or 7r are saddles is related to 
a fundamental property of the interacting solitons: while an effective mass, 
tur, of the two-soliton system corresponding to the radial degree of freedom 
R is positive, an effective mass of the phase degree of freedom is negative 
0- 
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Thus, these two types of the BS's are, respectively, unstable and stable, 
in the first approximation. In JOJ, it has also been demonstrated that the 
BS with ij; = 7r/2 is subject, in the next approximation, to an extremely 
weak instability, which transforms it into a very slowly unwinding spiral. 
However, it was also demonstrated that, even if this next-order instability 
can be observed, it does not destroy the BS, but, instead, makes it dynamical, 
with R and ip very slowly oscillating in a limited range. Note that the same 
mechanism gives rise, in the ID case, to (almost) stable chains of the bound 
solitons; in the 2D and 3D cases, a new possible pattern is a lattice of the 
bound solitons. There may also exist "covalent soliton molecules", in the 
form of triangles and tetrahedrons in the 2D and 3D cases, respectively. 

In the absence of the dissipation, BS of mutually orbiting solitons is 
possible in the 2D and 3D cases (in the latter case, it is assumed that the two 
solitons move in one plane). Orbiting of incoherently interacting 2D solitons 
was experimentally observed in a photorefractive medium ||. Numerical 
simulations and analytical arguments presented in demonstrate that the 
orbiting BS states of the 2D solitons in the model are unstable. In 
the present class of the models, the orbiting BS cannot be stable either. 
Indeed, for the orbiting state the interaction potential fl29|) or ( |3T| ) must 
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be supplemented by the centrifugal energy E cl = ( y M 2 /2mR) R~ 2 , where M 
is the angular moment of the soliton pair, and mn is the above-mentioned 
effective mass. Thus, the net effective energy of the orbiting state is 

E cS = U D (R,^)+E cf = Cocostp iT (D ~ 1)/2 exp [-^/^%jR) + (m 2 /2m fl ) R~ 

(34) 

where, according to Eqs. fl29|) and ([H]), the constant Co depends on the 



dimension D and the soliton's spin s. It is easy to check that the effective 
energy ( |3~4"D gives rise to a stationary state with sin?/' = 0, Cocosip < 0, 
provided that M 2 is small enough. However, this stationary state always has 
d 2 E c g/dR 2 < 0, i.e., it is a maximum of the effective energy, consequently, 
the orbiting BS is unstable against variation of R. Moreover, one can check 
that the same state always has d 2 E eS /dip 2 > O.With regard to the above- 
mentioned m^p < 0, this BS is also unstable against variation of ip. 

In conclusion, a general method to find the effective potential of interac- 
tion between two-dimensional and three-dimensional solitons was elaborated, 
including the case of the two-dimensional vortex (spinning) solitons. The 
method is based on calculation of the overlapping term in the full Hamil- 
tonian of the system. The main technical point that makes the calculation 
possible is that the bulk integral reduces to a surface one, and, in the lowest- 
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order approximation, the final expression does not contain the auxiliary ra- 
dius of the overlapping region. The result applies to spatial solitons and 
"light bullets" (spatiotemporal solitons) in nonlinear optics (in the model 
with the quadratic nonlinearity, the interaction between the "bullets" may 
be spatiotemporally anisotropic). The interaction potential predicts that an 
orbiting bound state of two solitons exists, but is always unstable. In the 
presence of weak dissipation and gain, the effective potential can also be de- 
rived, giving rise to bound states of the solitons (both unstable and almost 
stable) similar to those recently studied in the one- dimensional model. 
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5 FIGURE CAPTION 

Fig. I. The two-soliton configuration in the two- and three-dimensional 
cases (in the latter case, the figure shows the central cross section of the 3D 
configuration). The points 1 and 2 are centers of the solitons. 
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Fig. 1 



